Note: This example does not yield the exact same results as https://niivue.com/demos/features/denoise.html. This notebook does, however, show how one can change the volume img data and save it.

from pathlib import Path

from ipyniivue import download_dataset

BASE_API_URL = "https://niivue.com/demos/images/"
DATA_FOLDER = Path("images")

# Download data for example
download_dataset(
    BASE_API_URL,
    DATA_FOLDER,
    files=[
        "otsu.nii.gz",
    ],
)

otsu.nii.gz already exists.
Dataset downloaded successfully to images.

try:
    from scipy.ndimage import grey_dilation, median_filter
except ModuleNotFoundError:
    !pip install scipy
    from scipy.ndimage import grey_dilation, median_filter

import ipywidgets as widgets
import numpy as np

from ipyniivue import NiiVue, Volume
from ipyniivue.constants import DragMode, ShowRender
from ipyniivue.utils import find_otsu

# Setup niivue object
nv = NiiVue(
    back_color=(0.9, 0.9, 1, 1),
    show_3d_crosshair=True,
)

nv.load_volumes([Volume(path=DATA_FOLDER / "otsu.nii.gz")])

imgRaw = None


@nv.on_image_loaded
def on_image_loaded(volume):
    """Handle image loaded."""
    global imgRaw
    imgRaw = volume.img.copy()
    print("Image loaded. Original image data stored.")


@nv.on_location_change
def handle_intensity_change(location):
    """Handle location change."""
    with intensity_output:
        intensity_output.clear_output()
        print(location["string"])


nv.opts.multiplanar_show_render = ShowRender.ALWAYS
nv.opts.yoke_3d_to_2d_zoom = True
nv.set_interpolation(True)
nv.set_clip_plane(0.2, 0, 120)

# Create UI elements

denoise_button = widgets.Button(description="Denoise")
save_button = widgets.Button(description="Save")
cancel_button = widgets.Button(description="Cancel")
apply_button = widgets.Button(description="Apply")

otsu_select = widgets.Dropdown(
    options=[
        ("Very Heavy 3:4", 1),
        ("Heavy 2:3", 2),
        ("Medium 1:2", 3),
        ("Light 1:3", 4),
        ("Very Light 1:4", 5),
        ("None", 6),
    ],
    value=5,
    description="Reduce dark noise:",
)

dilate_check = widgets.Checkbox(value=True, description="Dilate dark")
denoise_check = widgets.Checkbox(value=False, description="Denoise")
dark_check = widgets.Checkbox(value=True, description="Clip Dark")
intensity_output = widgets.Output()

# Functions to handle UI events


def on_dark_check_change(change):
    """Handle dark check change."""
    nv.opts.is_alpha_clip_dark = change["new"]


on_dark_check_change({"new": dark_check.value})


def on_drag_mode_change(change):
    """Handle drag mode change."""
    if change.new == "none":
        nv.opts.drag_mode = DragMode.NONE
    elif change.new == "contrast":
        nv.opts.drag_mode = DragMode.CONTRAST
    elif change.new == "measurement":
        nv.opts.drag_mode = DragMode.MEASUREMENT
    elif change.new == "pan":
        nv.opts.drag_mode = DragMode.PAN


drag_mode_dropdown = widgets.Dropdown(
    options=[
        ("Drag pan/zoom", "pan"),
        ("Drag contrast", "contrast"),
        ("Drag measurement", "measurement"),
        ("Drag none", "none"),
    ],
    value="pan",
    description="Drag Mode:",
)
drag_mode_dropdown.observe(on_drag_mode_change, names="value")


def on_denoise_button_clicked(b):
    """Denoise."""
    if denoise_options.layout.display == "none":
        denoise_options.layout.display = ""
    else:
        denoise_options.layout.display = "none"


denoise_button.on_click(on_denoise_button_clicked)


def on_cancel_clicked(b):
    """Cancel/reset operation."""
    global imgRaw
    nv.volumes[0].img = imgRaw.copy()


def on_save_clicked(b):
    """Save image."""
    nv.volumes[0].save_to_disk("denoised.nii")


# Process img functions


def denoise_image(img_raw, nx, ny, nz):
    """Denoise image."""
    img_raw = img_raw.reshape(nz, ny, nx)
    img_smoothed = median_filter(img_raw, size=3)
    return img_smoothed.flatten()


def dilate_image(img_raw, nx, ny, nz):
    """Dilate image."""
    img_raw = img_raw.reshape(nz, ny, nx)
    img_dilated = grey_dilation(img_raw, size=(3, 3, 3))
    return img_dilated.flatten()


def process_volume(change=None):
    """Process volume."""
    global imgRaw
    if imgRaw is None:
        print("Image not loaded yet.")
        return

    level = otsu_select.value

    # Reload original image
    img_raw = imgRaw.copy()
    nx = int(nv.volumes[0].dims[1])
    ny = int(nv.volumes[0].dims[2])
    nz = int(nv.volumes[0].dims[3])

    if level in [5, 1]:
        otsu = 4
    elif level in [4, 2]:
        otsu = 3
    else:
        otsu = 2

    thresholds = find_otsu(nv.volumes[0], mlevel=otsu)
    if len(thresholds) < 3:
        print("Threshold calculation failed.")
        return

    threshold = thresholds[0]
    if level == 1:
        threshold = thresholds[2]
    elif level == 2:
        threshold = thresholds[1]

    mn = np.min(img_raw)
    if level > 5:
        threshold = mn

    if denoise_check.value:
        img_raw = denoise_image(img_raw, nx, ny, nz)

    if dilate_check.value and level < 6:
        imgMx = dilate_image(img_raw, nx, ny, nz)
    else:
        imgMx = img_raw.copy()

    # Apply threshold
    mask = imgMx < threshold
    img_processed = img_raw.copy()
    img_processed[mask] = mn

    img_processed = img_processed.astype(imgRaw.dtype)

    nv.volumes[0].img = img_processed


# Attach handlers

apply_button.on_click(process_volume)
cancel_button.on_click(on_cancel_clicked)
save_button.on_click(on_save_clicked)
dark_check.observe(on_dark_check_change, names="value")
otsu_select.observe(process_volume, names="value")
dilate_check.observe(process_volume, names="value")
denoise_check.observe(process_volume, names="value")

# Display all

denoise_options = widgets.VBox(
    [
        otsu_select,
        dilate_check,
        denoise_check,
        apply_button,
    ]
)
denoise_options.layout.display = "none"

ui_controls = widgets.VBox(
    [
        widgets.HBox([denoise_button, save_button, cancel_button]),
        denoise_options,
        dark_check,
        drag_mode_dropdown,
    ]
)

display(ui_controls)
display(nv)
display(intensity_output)